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ABSTRACT 


A major step in a most probable point (MPP)-based method for 
reliability analysis is to determine the MPP. This is usually 
accomplished by using an optimization search algorithm. The minimum 
distance associated with the MPP provides a measurement of safety 
probability, which can be obtained by approximate probability 
integration methods such as FORM or SORM. The reliability sensitivity 
equations are derived first in this paper, based on the derivatives of the 
optimal solution. Examples are provided later to demonstrate the use of 
these derivatives for better reliability analysis and reliability-based 
design optimization (RBDO). 
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I. INTRODUCTION 


Generally, sampling and approximate integration are two commonly used approaches for 
reliability analysis of a general response or a limit-state equation. The Monte Carlo simulation 
plays a key role in the sampling approach, whereas the concept of MPP does the same in the 
approximate integration approach. Derivatives of quantities resulting from reliability analysis 
can also be computed by either the sampling approach [1] or the approximate integration [2]. The 
focus of this paper is on the sensitivity analysis pertaining to the approximate integration 
approach. 

As revealed by the Taylor’s series expansion, derivatives of a function provide a means 
to approximate the function beyond the current design point. Therefore, such derivatives can be 
used to answer the ‘what-if’ questions in a design environment. This important feature of the 
derivatives leads to the development of a discipline, called sensitivity analysis, which seeks 
efficient computational methods to compute the derivatives. Many works in reliability 
engineering, particularly in reliability-based design optimization, have given special attention to 
study the related derivatives [1-9]. However, some of these derivatives [2] are computed 
according to the functional relationship defined by the limit-state equation. Such derivatives are 
called behavior derivatives [10,11]. Since the MPP can be obtained as a result of an optimization 
process, the probabilistic derivatives can also be viewed as the derivatives of the optimal solution 
or optimum sensitivity derivatives [10,11]. In other words, the computation of such derivatives 
should involve not only the function of the limit- state equation but also the Kuhn-Tucker 
Necessary Conditions at the MPP. 

The main goal of this paper is thus to present a procedure that computes the probabilistic 
derivatives as derivatives of an optimal solution. Furthermore, examples are used to demonstrate 
the applications of such optimum sensitivity derivatives to form better procedures for reliability 
analysis and reliability-based design optimization. 

The paper is organized into three main parts. The first part, Sections II and III, reviews 
the procedures to derive optimum sensitivity derivatives and to perform the reliability analysis. 
The second part, Sections IV and V, derives the sensitivities of the results of reliability analysis 
with which a new MPP-based reliability analysis method is established. Numerical verifications 
are given in Section VI where simple academic examples are used to demonstrate the use of the 
newly devised method for reliability analysis and design optimization. 

II. DERIVATIVES OF CONSTRAINED OPTIMAL SOLUTIONS 

A typical engineering optimization problem seeks the solution, which achieves the best 
performance among a given set of problem parameters within the limits of the available 
resources. The best performance is usually measured by a targeted objective, while the limited 
resources are usually prescribed by constraints. Mathematically, an optimization problem can be 
symbolically represented as [12] 
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min f(b, p) 

beR n 

subject to 

h i {b,p) = 0, i = 1,2, . . . m /? 
8j{b,p)< 0, j = l,2,---m g 


(P-1) 


where b e R” is the design variable vector which is unknown, p is the given problem parameter, 
Jlb.p) is the objective and hj(b,p ) and gj(b,p ) are the equality and inequality constraints, 
respectively. The integers, /«/, and m g , represent the numbers of equality and inequality 
constraints, respectively. Note that the problem parameter can be a vector; however, for 
simplicity, it is taken as a scalar here. 

Let the optimal solution of Problem (P.l) be set as b e R 11 . It is then known [12] that, b *, 
satisfies the Kuhn-Tucker Necessary Conditions, 


^r(b\ p)=0 

db 


(1) 

h i (b*,p)= 0, 

i = 1,2, . . . 

(2) 

8j(b*,p)<0, 

j = l,2,...m g 

(3) 

A j8j(b\p)=0 , 

j = 1,2, . . . m g 

(4) 

^>0, 

j = 1 , 2 ,... m g 

(5) 


where L is the Lagrangian, defined as 

m h m g 

L^f + Y j r i h i +Y J A j g j ( 6 ) 

‘ j 

and the Lagrange multipliers, r, and A, are uniquely defined if the gradients of the objective and 
the constraints are linearly independent at b . According to Eqs. (3-5), the inequality constraints 
at the optimal solution, b , can be classified into two possibilities; tight or non-tight. The tight 
constraints follow the relations 


8j(b*,p)=0, j e m g 

Aj >0, j e in g 


(7) 

( 8 ) 


while the non-tight ones follow a different set of relations 
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(9) 

( 10 ) 


gi(b*,p)<0, j£m g 

Aj= 0, j<£m g 


where m is a collection of tight constraints. Consequently, the Lagrangian of Eq. (6) can be 

o ' 

simplified as 


m h 

L =f+Tj r i h i + (ii) 

i jem g 


The solution of an optimization problem should then include b as well as the Lagrange 
multipliers as unknowns. Note that the Lagrange multipliers associated with tight inequality 
constraints are subjected to sign restriction, as indicated by Eq. (8), whereas those associated 
with equality constraints are not. The Kuhn-Tucker Necessary Conditions clearly show that the 
optimal solutions, b and the Lagrange multipliers, will be changed if p is changed. Thus, the 
optimal solutions are functions of problem parameter, p, and their derivatives with respect to p 
are called optimum sensitivity derivatives. 

II.l Optimum Sensitivity Derivative 


An optimum sensitivity derivative is the limit of the difference between two neighboring 
optimal solutions that satisfy respective Kuhn-Tucker Necessary Conditions with different p’s. 
Therefore, the optimum sensitivity derivatives can be obtained by differentiating the Kuhn- 
Tucker Necessary Conditions, Eqs. (1-2) and Eqs. (7, 11) with respect to p as 


d 2 L db* 
db 2 dp 


IXdh dr, ^dgjdAj d 2 L 
i db dp db dp dbdp 


( 

v 

r 

V 


dh; 

~db 


w T 

db 


db * dh 
~dt> + 3p = °’ 


dp dp 


i = 1,2, . . . m h 


j e m 
J 8 


( 12 ) 

(13) 

(14) 


where the Hessian of L and 


d 2 L 

dbdp 


are given as 


d 2 L d 2 f 44 d 2 h; ^ . d 2 gj 
db 2 db 2 + ^ n db 2 + % Xj db 2 


(15) 
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( 16 ) 


d 2 L d 2 f A d 2 h t ^ 5 d 2 g 
dbdp dbdp j dbdp J dbdp 


For clarification, Eqs. (12-14) can be put in a matrix form as 


d 2 L 

dh 


C -i. > 

db 


d 2 L 

db 2 


db 

db 

dp 


dbdp 

fdhj 

T 

0 

0 

dr 


dh 




> — — < 

— 

U* J 




dp 


dp 

ra^'1 

T 

0 

0 

dA 


¥ 

UJ 


dp 


dp 


(17) 


where the summation notations in Eq. (17) are dropped for simplicity. Note that such 
differentiability is valid only when the tight constraint set remains the same in the neighborhood 
of p. In other words, Eq. (17) is valid only if the equations and the numbers of the tight 
constraints and their associated Lagrange multipliers as prescribed by Eqs. (7-8) remain 

s|s _ 

unchanged around p. Otherwise, the optimal solutions b and A are not differentiable with respect 
to the parameter at the given value of p. 

The dimension of the leading coefficient matrix in Eq. (17) is m = m h + m g and the 

unknowns are db* / dp , dr I dp and dA/dp for a single parameter, p. Though Eq. (17) has a 
unique solution as the leading coefficient matrix is non-singular, it is difficult to construct for 
general engineering applications. This is because its leading coefficient matrix involves second- 
order derivatives of L, such as Eqs. (15-16). Nevertheless, one particular type of the optimum 
sensitivity derivative, the derivatives of the objective, is much easier to obtain than the unknowns 
of Eq. (17). The derivation of such derivatives follows. 

II. 2 Optimum Sensitivity Derivatives of the Objective 

>(: 

Differentiation of the objective at b with respect to p yields 


df_ 

dp 


df 

= — + 
dp 


¥ 

V db J 


db * 
dp 


The derivative, df/db , in the above equation can be replaced by the derivatives of the 
constraints, as stated by Eqs. (1, 11), 
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(18) 


d[_ 

dp 


dp 



dp 


-14 


jem g 


\ M J 


db* 

dp 


The terms associated with db / dp in the right-hand side of Eq. (18) can be replaced by the 
derivatives of the tight constraints. That is, according to Eqs. (13-14) and Eq. (11), the end result 
is given as 


df df 44 dh ^ , dgj dL 

dp dp f dp je „, s dp dp 


(19) 


Equation (19) states that the derivative of the optimum objective with respect to a parameter is 
equal to the derivative of the Lagrangian with respect to the same parameter at the optimal 
solution. The equation gives a first order approximation of the change in the optimum objective 
due to the change in the problem parameter. 

III. RELIABILITY ANALYSIS 

Given a response condition of random variables, X, reliability analysis is interested in 
finding the probability of failure of such a condition. For example, a constraint is defined by the 
response condition 


G(A)<0 


( 20 ) 


which gives the probability of failure as 


P f =P{G{x)> 0 ) ( 21 ) 

The reliability is then given by 

R = \-P f (22) 

The limit-state equation is defined by the failure surface, which is given by 

G(x) = 0. (23) 

Note that in this study, the random variables, X, are assumed to be statistically independent and 
normally distributed. If they are not, one should take proper procedures to find their suitable 
equivalences [2] before reliability analysis. 

One approach to solve Eq. (21) is the approximate probability integration method whose 
goal is to find a rotationally invariant reliability measurement, with which First-Order Reliability 
Method (FORM) or Second-Order Reliability Method (SORM) is developed. Most of these 
methods may be classified into two groups; the Reliability Index Approach (RIA) and the 
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Performance Measure Approach (PMA) [6, 7]. As an example, the HL-RF algorithm [13-16] is a 
method of RIA, whereas the AMV+ [1, 17] is one of PMA. 

In RIA, the objective is to find the first-order safety reliability index, /?, which is equal to 
the shortest distance between the origin of the C-space and the failure surface. The contact point 
on the failure surface is called the Most Probable Point (MPP). The u is the reduced variable 
vector, whose component is defined by 


where jUj and cr l are the mean and the standard deviation of the corresponding X,. The limit-state 
equation is thus rewritten as G(u, <7, ju) = 0. Mathematically, the reliability index can be viewed 
as the objective of an optimization problem with the limit-state equation as the equality 
constraint. 


min if 

subject to «eR"" " 
G(u) = 0 


(P.2) 


where the reduced variable vector, u, is considered as the design variable vector. At its optimal 
solution, u* , one finds the reliability index, /? = ||m || , which yields the first-order probability of 

failure as Pj = c P(-/3). 

In PMA, the objective is to compute the first-order probabilistic performance measure, 
G * . It is defined as the offset of the performance, G(X), so that the shortest distance between the 

limit-state equation, G(X) - G p = 0 and the origin of the U - space is equal to a given target 

reliability service, <f>(- jB o ) . G p can also be found as the smallest value of G that is tangent to 

the targeted reliability surface, represented by a sphere constraint, ||w|| = J3 n . Mathematically, the 

first-order probabilistic performance measure is obtained as the objective of an optimization 
problem as 

minG(i<) ] 

subject to “**" >• (P-3) 

F = A, ! 


where the constraint requires the solution to achieve the targeted reliability index. 

The optimization formulation of the RIA is very similar to that of the PMA. In fact, at 
their optimal solutions, their Kuhn-Tucker Necessary Conditions yield similar relations; 
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for the RIA, 


(24) 


u 
I u 


+ A 


dG 

du 


= 0 


and, for the PMA, 


i u dG _ 


ll 


du 


(25) 


It is noted that their Lagrange multipliers can be either positive or negative, since they 
correspond to equality constraints. Furthermore, their Lagrange multipliers exhibit a reciprocal 

>j: :js 

relationship, if b is the same for both problems. However, their optimal solutions, b , are usually 
not the same, unless the targeted reliability index in PMA is exactly the same as the reliability 
index obtained in RIA. Since Hi/ II = J3* in RIA and ll/i II = j3„ in PMA, Eqs. (24-25) give 
convenient means to compute the respective Lagrange multipliers as 


A 




dG 

\du j 


u 


and 


dG 

; _ V du j 

p ~ aT 


(26) 


(27) 


IV. SENSITIVITIES OF RELIABILITY ANALYSIS RESULTS 

Since the standard reliability methods, RIA or PMA, are all formulated as constrained 
optimization problems, Eqs. (17) and (19) can be directly applied here to calculate their optimum 
sensitivity derivatives. Note that in Problems (P.2) and (P.3), the reduced variables, u, are treated 
as the design variables, whereas the standard deviations a, the mean values jll of the random 
variables or the targeted f3o in PMA are treated as the problem parameters. Furthermore, the 
optimum sensitivity derivatives of Problems (P.2) and (P.3) are readily available, as they involve 
only equality constraints. 

IV. 1 Problem (P.2) for RIA 

The sensitivity of the reliability index, p , with respect to a problem parameter can be 
obtained from Eq. (19) as 


dp * _ dp ' | dG dG 
dp dp ' dp ’ dp 


( 28 ) 
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where 90 /dp = 0 since 0 is a function of u only. To find the sensitivity derivatives of the 
optimum reduced variables, u*. an equation of the second order derivatives of 0 and the 
constraint, G(u) = 0, is then required, as indicated by Eq. (17), 


* * T 

u 

. d 2 G dG 
+ K t - r 

r * > 

du 


A 

r d 2 c Y 


du 2 du 

dp 

__ 

V 

0udp 

f 3G l 

T 

dA r 



dG 

v du j 

0 

dp 


_ 

dp 


(29) 


The random variable, X t , is related to its reduced form, «„ as X, = «, a, + //,. Thus, the 
components of the gradient and the Hessian of G(u ) or G(X) = 0 in Eq. (29) are given as 

dG _ f 5G ] _ \ dG ] 
du {dwj [dX, '] 

and 


d 2 G _ 

d 2 G 


d 2 G 

du 2 ' 

du i du j 


dx [ dx j 


( 7,(7 


IV. 2 Problem (P.3) ofPMA 

The objective function /and the equality constraint in this case are G(u, o. 0) and h{u,0) 
= (u T u) U2 - 0) = 0, respectively. The derivatives of the first-order probabilistic performance 
measure, G* p , with respect to the problem parameters are given by Eq. (19) as 


dG p _ dG p ^ ^ dh_ 

dp dp p dp 


(30) 


In particular, if the standard deviations and the mean values of random variables are considered 
as the problem parameters, then Eq. (30) is simplified as 



dp dp 


(31) 


as the constraint, h, is a function of u and 0 only. On the other hand, if the target 0 is selected 
as the problem parameter, the first term, dG * /dp , becomes zero and the second term, oh/ dp, 
equals -1. As a result, Eq. (30) is reduced to 
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( 32 ) 


dfi 0 



>(t 

The sensitivity derivatives of the optimum solution, u , of the PMA with respect to the problem 
parameters are given by Eq. (17) as 


du ' 


U^-i 

* 

u 

r * 'i 

du 

A, 

Jo 

dp 

u] 

T 

0 

d\ 

A) , 


dp 


d 2 G . d 2 h 

+ K 

du dp p du dp 
dll 

dp 


(33) 


where the relations, dh/du = u */J% and d 2 h/du 2 = I/J%, have been used in derivation. 

V. APPLICATIONS 


The optimum sensitivity derivatives derived above provide means to support the 
development of a better algorithm for reliability analysis as well as reliability-based design 
optimization (RBDO). 

V. 1 Gradients for RBDO 


A reliability-based design optimization usually involves random variables as design 
variables and reliability index or performance measure, depending upon its formulation, as its 
objective or constraints. The randomness of a design variable in a reliability-based design 
optimization is usually represented by its mean and standard deviation. Therefore, the mean and 
the standard deviation of the random variables can be directly modeled as the design variables. 
In this case, Eq. (28) or Eq. (31) provides necessary derivatives of the reliability index or the 
performance measure with respect to the mean or the standard deviation of a random variable to 
support any reliability-based design optimization algorithm. More specifically, in terms of the 
random variable, X,-, the derivatives of the reliability index of Eq. (28) with respect to the mean 
and the standard deviation are rewritten, respectively, as 


dp _ dG 
d J u i ~ r 3 A 


(34) 


and 


dp’ 

d<7; 


A 


dG 

ax, 


(35) 


On the other hand, the derivatives of performance measure of Eq. (31) with respect to the mean 
and the standard deviation are given as 
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( 36 ) 


and 


dG\ _ d G 
dll, dX, 

dG^ = dG_ * 
da, dX,. U ‘ 


(37) 


V.2 Search Algorithms for Reliability Analysis 


Equations (24-25) provide an important rule to guide the search of the reliability index or 
performance measurement that u is proportional to dG/du. Thus, the unit direction of u is given 
as 


dG/du 

dG/du\\ 


(38) 


In PMA, since the length of vector u is limited to be the given /3o, as prescribed in the 
constraint, a very simple recursive formula can then be devised to find the u as 


z+1 


A)«' 


(39) 


This algorithm is quite efficient. However, Choi and Youn experienced convergence difficulty of 
Eq. (39) when the limit-state equation exhibits concavity [7]. They then replaced the current unit 
direction, n\ in Eq. (39) by the average of the last three consecutive unit directions. Their 
investigation can be found in Reference 7. 

The search algorithm in RIA is more complicated than that in PMA, because 0 is not 
known in advance. Equation (39) is still valid; however, a correction step is needed to improve 
the estimate of j3. A recursive formula on u as well as j3 can then be set up as 

u M = j3 M n' (40) 

where 0 +l can be obtained through the Newton-Raphson’s iteration on the constraint, G (u ,+ ) = 
G(0 +t n ) = 0, as 


■AgY 2 

n 


v du J 


A p‘ =-G(p‘ii) 


(41) 
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Once P ' 1 is converged, u ,+l is then found by Eq. (40). One can then move to next iteration. This 
constitutes the basic algorithm of Hasofer-Lind method [13-14]. An improved version [2,15,16] 
of the Hasofer-Lind method is to reduce computational burden of the Newton-Raphson’s 
iteration by updating u and p in the same iteration as 

u i+1 =P‘n i +Ap‘n‘ (42) 


where P is obtained by 


P 


= u n 


and AP is obtained by solving Eq. (41) as 


A 0 = 


-G(p‘ri) 



These RIA methods, however, fail to converge, if the limit-state function is concave near 
the search point [7]. A more robust method for RIA is proposed hereafter, which is derived based 
on the PMA algorithm and its optimum sensitivity equation, Eq. (32). 


V.3 An Alternative Algorithm for RIA 


The key motivation of this new algorithm is the observation that the target Po of PMA is 
identical to the reliability index, 0, of RIA if the performance measure, G* , in PMA reaches 

zero value. To achieve a zero G * , the new algorithm uses Eq. (32) to estimate the amount of the 
change in /? 0 needed to reduce the non-zero value of G . Eq. (32) gives an estimate, 

AG* = -A p A p (43) 

Let AG = 0 -G* be the gap between zero and the current value of G* . The change in target Po 
that is required to achieve a zero G* is then estimated by Eq. (43) as 


A P 




(44) 


The updated reliability index, P + AP , will yield a new G p that is closer to zero, at least in the 
first order sense. Repeated use of Eq. (44) can lead the PMA search to G* = 0. In short, the new 
algorithm can be summarized as follows: 

Step 1: Start with an initial target p' = po and initial values of u = u () . 

Step 2: Follow the PMA procedure [7] to obtain the converged performance measure, G p . 

Step 3: Compute A P by Eq. (27), which is evaluated at the converged u and p l in Step 2. 
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Step 4: Update j3 ,+1 = /3' + A[3 ', where A(3 1 = G*‘ IA P . 

Step 5: Return to Step 2, until A(3 and G* achieve the required tolerances for convergence. 

The above algorithm is similar to the one presented by Du and Chen [18]. Both trace the 
MPP locus to locate the 0. Nevertheless, the former follows the tangential direction of the MPP 
locus, whereas the latter follows an extrapolated MPP locus. A MPP locus, represented by the 
curve of G* (/?), is shown in Fig. 1. It is also important to note that Eq. (44) is in the same form 

as Eq. (41). Nevertheless, Eq. (44) is only valid at the optimum solution of the PMA. Therefore, 
it represents a slope on the MPP locus, whereas Eq. (41) does not. 

Since j3 is restricted to be positive, the G* {0) curve does exhibit a discontinuity at /3 = 0, 

as shown in Fig. 1. Fortunately, the probabilities of failure for most of applications are pertaining 
to the tip of either the upper or the lower branch of the curve, so as to avoid the discontinuity. 

The detailed discussion of this matter is given in the next section. 

V.4 Infeasible Region: G(X) <0 or G(X) >0 

In a classical design optimization problem, a feasible design is usually defined by G(X) < 
0, which is associated with a failure probability of Pf(G(X) > 0). The PDF plot and the MPP 
locus of this case are shown in Fig. 2. However, in a reliability-based design optimization, a 
reliable constraint is usually defined by G(A) > 0 [2], which is associated with a failure 
probability of Pf( G(X) < 0). This case is shown in Fig. 3. Note that Pf(G(X ) < 0) occurs on the 
far left tail of the probability distribution of G(u) along the direction of u in the standard normal 
space, whereas Pf (G(X) > 0) on the far right tail. This inconsistency, also noted by 
Padmanabhan and Batill [19], can cause confusion in selecting a proper solution procedure for 
reliability analysis. For example, the reliability analysis formulation and procedure derived in 
Sections V.2 and V.3 are based upon the definition of failure probability, Pf (G(X) > 0). 
Modifications of these formulation and procedure for the case of Pf(G(X) < 0) are in order. 

Note that both cases, Pf(G(X) > 0) and Pf(G(X) < 0), enjoy the same limit-state equation, 
G(X) =0. Therefore, both failure probabilities can be calculated from the same reliability index. 
Since the reliability index is always positive as defined by 0 = Hull, the recursive equation of 
Eq. (39) or (40) are valid only when u and n are in the same direction, where n is defined by 
Eq. (38) as the steepest ascending direction at the given point on the G(u ) curve. In the case of Pf 
(G(X) > 0), u is indeed in the same direction of n, as shown in Fig. 4, whereas in the case of Pf 
( G(X ) < 0), u is not, as shown in Fig. 5. Thus, the basic equation, Eq. (39) or (40), is valid only 
for the case of Pf(G(X) > 0). Otherwise, a negative sign should be added to these equations. For 
example, Eq. (40) should be modified as 


u M =-0 +l n‘ (45) 

It then follows that A r and A P are positive for the case of Pf(G(X) < 0) and negative for the case 
of Pf(G(X) > 0). Consequently, the slope of the MPP locus, G* (/?), is negative for the former 
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case and positive for the latter, as indicated by Eq. (32). Figures 2.b and 3.b demonstrate the 
slopes of the MPP locus for both cases. 

>(t _ 

The difference in the slopes of G p ip) has an important implication in the formulation of 
reliability constraints. For example, if a reliability constraint requires the failure probability of 
constraint, G(X) > 0, to be less than a targeted value, Pq, the reliability constraint can be stated as 

P f (G(x)<0)<P 0 (46) 

The above constraint is equivalent to 


A s A (47) 

where the reliability index, p , is calculated by (p (-/>) = Pf(G(X) < 0) and the targeted 0 is 
calculated by <t> i-fio) = Pq. Eq. (47) implies, according to Fig. 3.b, 

o <g;(A). (48) 

On the other hand, if the reliability constraint requires 

P f (G(X) > 0) < P 0 . (49) 

The constraint is equivalent to, in terms of reliability indices, 

A></ (50) 

which is the same as Eq. (47), or, in terms of performance measurement, according to Fig. 2.b, 

g;m<o. ( 5i ) 

In the RIA-based RBDO, reliability constraints are expressed in the form of either 
Eq. (47) or (50), whereas in the PMA-based RBDO, reliability constraints are expressed in the 
form of either Eq. (48) or (51). 


VI. NUMERICAL EXAMPLES 

Numerical studies are conducted in this section to verify the equations derived in 
Sections IV and V and to demonstrate their applications to reliability analysis and reliability- 
based design optimization. The example limit-state equations used in this study are the ones 
studied in Ref. 7. One of them is a convex function, and the other is a concave one. Reference 7 
reported that the HL-RF method failed to converge in both examples and the AMV+ method 
performed very well in the convex example, but failed to converge in the concave one. A hybrid 
method was then proposed in Ref. 7 that is workable for both cases. If the limit-state equation is 
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convex, the hybrid method will be the same as the AMV+ method, which will take the gradient 
of G(X) at the current design point as the search direction. If the limit-state equation is concave, 
the hybrid method will take an average of the gradients evaluated at the last three consecutive 
design points as the new search direction. 

The performance of the new RIA developed in Section V 3 for reliability analysis is first 
examined in Section VI. 1. The new RIA first uses the hybrid method [7] to find the performance 
measure, G p , for the targeted /? 0 . If the value of G p does not converge to zero, the hybrid 

9{S _ 

method is then restarted to find the G p of the new p 0 updated by Eq. (44). The process is 
continued until G p reaches zero within a tolerance. Once convergence, the derivatives of the 
results of reliability analysis are readily available, according to Eqs. (28) and (31), to support 
RBDO. 


In Section VI. 2, two approaches are first considered for RBDO. One is a RIA-based 
RBDO which models its reliability constraints in the form of Eq. (47) or (50), the other is a 
PMA-based RBDO which models its reliability constraints in the form of Eq. (48) or (51). The 
new RIA is incorporated with a gradient-based optimization algorithm, called L inearization 
Method [20, 21] to form the RIA-based RBDO solution procedure. On the other hand, the hybrid 
method is incorporated with the Linearization Method to form the PMA-based RBDO. The 
performance of the RIA-based RBDO will then be compared with that of a PMA-based RBDO. 
It is found that the RIA-based RBDO takes more function and gradient evaluation for reliability 
analysis than the PMA-based RBDO does. However, the former is much more robust than the 
latter in the line search part of a design optimization procedure. This observation leads to the 
development of a new RBDO procedure that incorporates the concept of active reliability 
constraint to achieve better computational efficiency in RBDO. 

The suggested method will start with the hybrid method to calculated G p for the targeted 

_ >(: 

po for each of the reliability constraints. If G p is less than a given margin, the constraint is 
judged as active-in-reliability. That means, the current design point is near the boundary of the 
reliability constraint. The new RIA will only be called upon to calculate the reliability indices of 
those active-in-reliability constraints. Once the values of such constraints are computed, the 
active constraint set can then be determined to guide the search direction in an optimization 
algorithm. Note that two different active constraint sets are mentioned here. One is active in 
reliability analysis, and the other is active in design optimization. The former set is bigger than 
the latter. 


VI.l The New RIA and Optimum Sensitivities 

The new RIA is first applied to find the reliability index of the limit-state equation, 

G l {X l ,X 2 ) = -e [Xl ~ l) -X 2 + 10 (52) 

which is a convex function. The random variables, X\ and 20 are both normal and given as X\ ~ 
N ( 6.0, 0.8 ) and X 2 ~ N ( 6.0, 0.8 ). The targeted j3 0 is set to be 3.0 and the initial values of X\ 
and X 2 are assigned to be 6 to start reliability analysis. The new RIA is called converged, if the 
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changes of the performance measurement, G p * and the reliability index, p , are, respectively, less 
than 0.0005 and 0.001 between consecutive iterations. The results are summarized in Table 1. It 
is noted that the new RIA takes 7 iterations to find the performance measure, G p = -0.357, for a 
po of 3 and it takes additional 5 iterations to reach the converged reliability index, [3 =2.878. The 
probability of failure, Pf (Gj < 0), is then calculated as 0.201%, which is in a good agreement 
with 0.216%, the probability of failure calculated by Monte Carlo simulation with 200,000 
samples. Note that each of iterations reported in the table involves a function and a gradient 
evaluation. 


Table 1: Convergence History of the New RIA - Example 1 


Iteration No. 

Xi 

*2 

G! 

0 

A p 

1 

6.828 

8.252 

0.905 

3.0 

0 

2 

7.546 

7.835 

0.437 

3.0 

0.9835 

3 

8.076 

7.202 

-0.138 

3.0 

1.5019 

4 

8.271 

6.773 

-0.341 

3.0 

2.4332 

5 

8.310 

6.647 

-0.357 

3.0 

2.9593 

6 

8.317 

6.624 

-0.357 

3.0 

3.0734 

7 

8.318 

6.620 

-0.357 

3.0 

3.0918 

1 

8.228 

6.596 

-1.413-2 

2.884 

3.2189 

2 

8.214 

6.647 

-1.669-2 

2.884 

2.8480 

3 

8.211 

6.656 

-1.677-2 

2.884 

2.8110 

1 

8.206 

6.565 

-8.358-5 

2.878 

2.8104 

2 

8.206 

6.659 

-9.260-5 

2.878 

2.7917 


The second example limit- state equation studied here is a concave function, 

G 2 {X l ,X 2 )= [e (a8XH ' 2) + e (01x '--°- 6] - 5]/ 10 (53) 

where X\ and X 2 are normally distributed random variables and X\ ~ N ( 4.0, 0.8) and Xi~N (5.0, 
0.8 ). The targeted fi 0 is set to be 3.0. Again, the convergence tolerances are taken as 0.0005 and 
0.001, respectively, for G p and Aft. The results are summarized in Table 2. The new RIA takes 9 
iterations to reach the converged G p of 0.204 for a targeted fio value of 3. Additional 7 iterations 
are needed for the method to find the converged reliability index, 0 = 3.803, which gives a 
probability of failure, P f (G 2 < 0), a value of 7.16E-3%. The probability of failure calculated by 
Monte Carlo simulation is 5.5 E-3% with 2 million samples. 


16 



Table 2: Convergence History of the New RIA - Example 2 


Iteration No. 

A, 

X 2 

Gz 

P 

A p 

1 

2.989 

2.823 

0.225 

3.0 

0.000 

2 

2.348 

3.259 

0.234 

3.0 

0.290 

3 

3.073 

2.786 

0.238 

3.0 

0.305 

4 

2.539 

3.096 

0.209 

3.0 

0.286 

5 

2.710 

2.976 

0.204 

3.0 

0.302 

6 

2.561 

3.068 

0.207 

3.0 

0.299 

7 

2.765 

2.942 

0.205 

3.0 

0.302 

8 

2.677 

2.998 

0.204 

3.0 

0.297 

9 

2.690 

2.989 

0.204 

3.0 

0.300 

1 

2.368 

2.549 

2.714E-2 

3.681 

0.244 

2 

2.310 

2.588 

2.725E-2 

3.681 

0.244 

1 

2.341 

2.450 

9.676E-4 

3.803 

0.217 

2 

2.199 

2.548 

1.628E-3 

3.803 

0.212 

3 

2.412 

2.405 

3.003E-3 

3.803 

0.214 

4 

2.233 

2.523 

8.631E-4 

3.803 

0.210 

5 

2.287 

2.486 

4.513E-4 

3.803 

0.213 


To verify the optimum sensitivities of the performance measurement and the reliability 
index, given by Eqs. (31) and (28), respectively, central differencing is conducted for the above 
two examples. In both cases, the mean value of X\, JU\, and the standard deviation of A3, 0}, are 
considered as the design parameters. The perturbations of them are taken as 0.01 and 0.008, 
respectively, in finite differencing studies. The results summarized in Table 3 show that the 
analytical derivatives calculated by Eqs. (31) and (28) are in good agreement with those 
calculated by central finite differencing (C.D.). 

VI. 2 Reliability-Based Design Optimization 


A general reliability-based design optimization problem can be symbolically written as 


subject to 


min f{ji) 

fieR" 

Pf * P, 


(P.4) 
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Table 3a: Optimum Sensitivities of G 1 and G 2 with respect to jUi 


Case 

+A/U\ 

~A/U\ 

C. D. 

Analytical 

Gi p 

-0.3955 

-0.3208 

-3.7371 

-3.7369 

G 2p 

0.2059 

0.2018 

0.2058 

0.2072 

p: 

2.8663 

2.8902 

-1.1967 

-1.1964 

fi *2 

3.8095 

3.7956 

0.6928 

0.7042 


Table 3b: Optimum Sensitivities of G\ 

and G 2 with respect to a 2 


Case 

+A(7 2 

-Ao 2 

C. D. 

Analytical 

Gi p 

-0.3642 

-0.3518 

-0.7785 

-0.7762 

G 2p 

0.1976 

0.2101 

0.7834 

0.7825 

p; 

2.8759 

2.8806 

-0.2954 

-0.2954 

pi 

3.7771 

3.8284 

-3.2037 

-3.2270 


where the probability of failure, Pf, can be calculated by either Pf (G(X) < 0 ) or Pj (G(X) > 0 ), 
depending upon the definition of failure. 

The specific test problem studied here is also taken from Reference 7. The design 
variables are the means of two statistically independent and normally distributed variables; X\ ~ 
N (5.0, 0.8) and X 2 ~ N (6.0, 0.8), where the initial guesses of the design variables are set to be 
5.0 and 6.0, respectively. The objective and the constraints of the problem are respectively given 
as 


f(/i) =3 m 2 - 2 m H 2 + 3 n 2 ~ 


(54) 


and 


P f (C,(X )<())< P l0 ; i = 1, 2, 


(55) 
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where the limit-state equations, Gi and G?, are defined by Eqs. (52) and (53), and Pm and P 20 are 
the required probabilities. 

In an RIA-based RBDO, the statements of failure probability, Eq. (55), are rewritten in 
the form of reliability index as indicated by Eq. (47), 

fi i0 -fi*<0; i = 1, 2 (56) 

whereas, in a PMA-based RBDO, they are rewritten in the form of performance measurements as 
indicated by Eq. (48), 


-G p *(Pio) < 0; i = 1, 2 (57) 

The Linearization Method [20, 21] is used here as the optimization algorithm to search 
for the best design. The Linearization Method is a special method of Sequential Quadratic 
Programming Technique in which an approximate sub-problem is formed in each optimization 
iteration as 


min -Sb T Sb + l T Sb 

Sbe R" 2 

h, + 1 1 Sb = 0; i = 1 to m, 

subject to 

g j + lgj3b - 0, ye A 


(P.5) 


where Sb is the change in design and l B , l h , and l g are the gradients of the objective 

function, equality and inequality constraints, respectively. The constraint set, A, includes only 
those constraints that are e - active. That is, the constraint, gj, that satisfies the condition, 

gj + z< 0 (58) 

The sub-problem, (P.5), is a convex problem that can be solved by a dual method. As a 
result, the computational effort of the method depends upon the number of the active constraints 
rather than the number of design variables. The solution of the sub-problem (P.5), Sb, gives a 
favored direction to reduce the cost as well the constraint violations. The new design point is 
then sought to be 


b i+I =b l +aSb 

where a is the step size. The value of a is selected to guarantee that a reduction in the merit 
function, y/, is achieved by the new design. That is, 
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if/{b 1+1 ) = if/{b ' + aSb ) < y/{b ' ) 


(59) 


The merit function, y/. is defined as 


xj/ = f + NF 

where F is the maximal value of the violations and N is a positive value, selected to be greater 
than the sum of all the Lagrange multipliers, 


m h 

l+ZT 

i=l jeA 

where r, and A j are the Lagrange multipliers corresponding to the equality and inequality 
constraints, respectively. Thus, Eq. (59) ensures that the new design will achieve a reduction in 
the cost as well as the constraint violations. It has been proven [21] that the Linearization 
Method reaches a local minimum if the Li-norm of the search direction approaches to zero, i.e., 

m<s o ( 60 ) 

The above equation serves conveniently as a convergence criterion. 

The convergence histories of the PMA-based and the RIA-based RBDO are tabulated in 
Tables 4 and 5. Both methods drove the L 2 -norm of the search direction to an acceptable level at 
the end of the optimization and they arrived at similar local minima. As expected, the PMA 
method is more efficient than the new RIA method for reliability analysis. The PMA-based 
method takes 131 and 39 function and gradient evaluations in total and the RIA-based method, 
155 and 111, respectively, for reliability analysis of limit-state equations, G i and G 2 . However, 
the PMA-based method takes 336 and 96 function and gradient evaluations for the purpose of 
line search, whereas the RIA-based method takes only 6 and 8. The line search factor makes the 
RIA-based method a better one for overall performance for this particular RBDO example. 

Note that each line search step in an optimization iteration requires a new reliability 
analysis. In the current study, the converged random variables, u , of the current line search step 
are kept as the initial values to start the next line search step. Furthermore, the final [3 and u at 
the end of the current optimization iteration are also kept as the initial values to start the next 
optimization iteration. 

The traditional design optimization procedure does not work well for RBDO. This is 
because the traditional procedure requires completed reliability analyses before judgment can be 
made regarding the activeness of the constraints. For example, the second constraint in the last 
example is never involved in the design optimization process. However, the RIA-based RBDO 
still performs 119 function and gradient evaluations to eventually determine the second 
reliability constraint is indeed not active. To reduce such unnecessary computational burden, the 


20 



Table 4. Results of PMA-based RBDO 


Iteration No. 

Cost 

Gpi 

Gp, 

Mi 

Mi 

No. of 
Analyses 
(Gpi/Gpi) 

IIAII 

Step-Size 

a 

No. of 
Analyses for 
Line Search 
{GpxtGpi) 

1 

123. 

0.971 

1.429 

5.000 

6.000 

9/4 

0.2634 

1.0 


2 

114.81 

0.779 

1.676 

4.850 

5.783 

9/3 

0.2545 

1.0 


3 

107.17 

0.618 

1.908 

4.704 

5.575 

9/3 

0.2457 

1.0 


4 

100.00 

0.481 

2.126 

4.562 

5.375 

9/3 

0.2373 

1.0 


5 

93.40 

0.366 

2.333 

4.423 

5.182 

9/2 

0.2292 

1.0 


6 

87.21 

0.267 

2.530 

4.288 

4.997 

9/2 

0.2213 

1.0 


7 

81.43 

0.182 

2.178 

4.157 

4.828 

7/2 

0.2133 

1.0 


8 

76.04 

0.110 

2.899 

4.0303 

4.647 

7/2 

0.2066 

1.0 


9 

71.00 

0.0469 

3.071 

3.906 

4.481 

7/2 

0.1995 

0.5 

7/2 

10 

68.63 

0.0188 

3.153 

3.845 

4.402 

7/2 

0.1962 

0.25 

7/2+7/2 

11 

67.48 

0.00556 

3.193 

3.816 

4.363 

7/2 

0.1948 

0.0625 

(7/2)*4 

12 

67.20 

0.00232 

3.203 

3.808 

4.348 

7/2 

0.1940 

0.03125 

(7/2)*5 

13 

67.06 

0.00071 

3.208 

3.804 

4.348 

7/2 

0.1939 

0.0078125 

(7/2)*7 

14 

67.03 

0.00031 

3.210 

3.804 

4.347 

7/2 

0.1938 

0.0039062 

(7/2)* 8 

15 

67.00 

0.00011 

3.210 

3.803 

4.347 

7/2 

0.1938 

0.0019531 

(7/2)*9 

16 

67.00 

1.38-06 

3.211 

3.803 

4.346 

7/2 

0.0059 

0.0002441 

(7/2)* 12 

17* 

67.00 

1.38-06 

3.211 

3.803 

4.346 

7/2 

0.0059 

- 



At the 17 th iteration, Sb\ = 0.0045, and Sb 2 = -0.00374 



Table 5. Results of RIA-based RBDO 


Iteration 

No. 

Cost 

A 

A 

Mi 

Mi 

No. of 
Analysis 

(A/A) 

IIAII 

Step-Size 

a 

No. of 
Analysis 
for Line 
Search 

(A/A) 

1 

123. 

5.549 

4.070 

5. 

6. 

45/22 

0.2635 

1.0 


2 

114.81 

5.221 

4.332 

4.850 

5.783 

22/9 

0.2545 

1.0 


3 

107.17 

4.904 

4.583 

4.704 

5.575 

26/8 

0.2457 

1.0 


4 

100.00 

4.595 

4.823 

4.562 

5.375 

15/8 

0.2373 

1.0 


5 

93.41 

4.298 

5.053 

4.423 

5.182 

15/8 

0.2292 

1.0 


6 

87.21 

4.011 

5.274 

4.288 

4.997 

9/8 

0.2214 

1.0 


7 

81.43 

3.738 

5.486 

4.157 

4.818 

20/8 

0.2138 

1.0 


8 

76.04 

3.473 

5.691 

4.030 

4.647 

21/8 

0.2066 

1.0 


9 

71.00 

3.209 

5.888 

3.906 

4.481 

12/8 

0.1995 

0.5 

6/8 

10 

68.63 

3.088 

5.982 

3.845 

4.401 

9/8 

0.0699 

1.0 


11 

67.00 

2.995 

6.050 

3.802 

4.347 

4/8 

0.0051 

1.0 


12 

67.08 

3.002 

6.050 

3.801 

4.352 

4/4 

0.0062 

1.0 


13* 

67.04 

3.002 

6.056 

3.795 

4.354 

1/4 

0.0049 

— 



*At the 13 th iteration, Sbi = 0.0025, Sb 2 = -.0.00425 


concept of active-in-reliability is implemented here. If a reliability constraint, defined by Eq. (57) 
satisfies the following condition, 

G' p (A)<e c (61) 

it is called active-in-reliability, because the design point is close to the boundary of the reliability 
constraint. The value of e G in Eq. (61) is small and positive. A small change in design can move 

the constraint from feasible to infeasible. The constraint of Eq. (61) can be monitored efficiently 
with the hybrid method, a PMA, which is the first step of the new RIA. Once the constraint is 
determined to be active-in-reliability, the new RIA can then be continued to finalize the 
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constraint value in the form of Eq. (56). With this new implementation, the last example is 
repeated here and the results are summarized in Table 6. The savings of computation on the 
second constraint is quite evident. In this example, the value of e G is specified to be 0.5. 

Table 6. Results of PMA/RIA-based RBDO 

No. of 


Analysis 

No. of for L ine 

Iteration Analyses Step-Size Search 


No. 

Cost 

A 

G p 2 

Mi 

Mi 

(G* pl /G p2 ) 

ll&ll 

a 

(G* pl /G p2 ) 

1 

123.00 

0.971 

1.429 

5.000 

6.000 

9/4 

0.2635 

1.0 


2 

114.81 

0.779 

1.676 

4.850 

5.783 

2/3 

0.2545 

1.0 


3 

107.17 

0.618 

1.908 

4.704 

5.575 

2/3 

0.2457 

1.0 


4 

100.05 

4.596 

2.127 

2.562 

5.375 

22/3 

0.2373 

1.0 


5 

93.40 

4.298 

2.333 

4.423 

5.182 

17/2 

0.2292 

1.0 


6 

87.21 

4.016 

2.531 

4.288 

4.997 

22/2 

0.2214 

1.0 


7 

81.43 

3.750 

2.719 

4.157 

4.818 

9/2 

0.2138 

1.0 


8 

76.04 

3.467 

2.899 

4.030 

4.647 

6/2 

0.2066 

1.0 


9 

71.00 

3.210 

3.071 

3.906 

4.481 

6/2 

0.1995 

1.0 


10 

68.64 

3.087 

3.153 

3.845 

4.402 

9/2 

0.0699 

0.5 

9/2 

11 

67.00 

2.996 

3.210 

3.802 

4.347 

4/2 

0.00501 

1.0 


12 

67.08 

3.002 

3.205 

3.801 

4.352 

4/2 

0.006109 

1.0 


13** 

67.04 

3.002 

3.203 

3.796 

4.354 

1/2 

0.004803 

.. 



* Starting from the 4 th iteration, the value is J3\ instead of Gp \ . 

**At the 13 th iteration, Sb\ = 0.00243, and = -0.00414 
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